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The purpose of the grant was to provide support to the author so that he 
could collaborate with Professor Gene H. Golub at Stanford University during 
the sunmer of 1978 on numerical problems relating to the Drazln inverse and 
applications thereof. In particular* the objective was to try to develop a 
useful numerical algorithm for the Drazln inverse and to analyze the numerical 
aspects of the applications of the Drazln inverse relating to the study of 
homogeneous Markov chains anJ systems of linear differential equations with 
singular coefficient matrices. 

It is felt that all objectives were accomplished with a measurable degree 
of success. 

A Stable Algorithm For the Drazln Inverse ; During the author’s stay at Stanford, 
the following algorithm was derived and reported on. The algorithm seems to be 
the best (in the sense of stability) in current existence. The algorithm is 
based on the work of Golub and Wilkinson (Ill-Conditioned Elgensystems and the 
Computation of the Jordan Form, SIAM Rev. , Vol. 18, Oct. 1976, pp. 578-619.) and 
is described below. 

The first half of the algorithm determines the index of a square matrix 
and returns the "index decomposition" of the matrix. 


I. The Index Decomposition. 

(0) Given A _ 0, apply the singular value decomposition to obtain 


A - U 


E 0 

rxr 


V. If r ■ n, then Ind(A) * 0 (l.e., A is nonsingular) 


D -1 * _i * 

and A - A ■ V E U and process stops. If r < n, then continue. 
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(1) Form VAV - VU 


Z 0 
0 0 


A. 


^(1) 

( 1 ) 

21 


. If AjJ^ - 0, then fP - 0. 


If a||^ ^ 0, apply the singular value decomposition to to obtain 


(1) ^1 ^ 

op 


where A^P Is n, x n, and E, Is r, x r, . Define P, “ V so that 


*11 


PlAPi 


‘1 


1 ‘1‘ 


^(1) 

^11 

,(1) 

21 


If r^ ■ then Ind(A) - 1 and Go To II. If < n^, then continue. 

.( 2 ) . 


(2) Form V 


A<“v * - V U ° 

1 11 1 11 [o 0 _ 


‘11 

,(2) 

^21 


If - 0, then 


( 2 ) 


A° ■ 0. If ^ 0, apply the singular value decomposition to 


to obtain 


( 2 ) 


11 




where 


( 2 ) 


Aj^' is ti2 X n^ and E^ is V2 ^ ^2* Define 


so that ^2^2 form 


P2AP2 


( 2 ) 

11_ 

( 2 ) 

21 

( 2 ) 

31 


A^2) 0 

^32 " 


If V 2 • n 2 » then Ind(A) » 2 and Go To II. If r 2 < n 2 » then continue. 
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( 3 ) ****** 


(k) 


fk-1) * 

Vi 


(k). 


(k) 

11 

(k) 

21 


. If ■ 0, then - 0. If 


(k) 


A^^V 0» apply singular value decomposition to A^^ to obtain 


o] 


(k) 

where A^^ is x and is r^ x Tj^. Define 


X-i “Ip 

.0 ij 


so that now has the form 



A(k) 

^11 

i” 

9 . . . . 

• 0 


A(k) 

^21 

1 

1 0 
1 

0 • • • • 

• 0 

■ 

,(k) 

31 

i 

j 32 

0 • • • • 

• 0 


A(k) 

41 

1 *42 

aC') . 
43 



,(k) 

\+l,l 

1 A^k) 

1 ^+1,2 

A.(k) . 

\+l,3 

• 

• 0 

■ n^, then Ind(A) ■ k and 

GO TO II. 

» 'k 


This process can continue for only a finite number of steps. When the process 

* 

terminates, the final value of is called the index decomposition of A. 
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* * 

If PAP Is the index deconposltlon, then F Is unitary and PAP has the 


form 



where N la lower triangular with O's on the diagonal (and hence nilpotent) and 
B la either nonaingular or a zero matrix. If B ° 0, then *0. If 6 is 
nonsingular, then 

D *r 

A*^ - P P 

.X 0. 

where X is computed by part II of the algorithm. 


II. Computation of the X Factor. 

Partition X and -CB ^ by rows so that 


' 


“ 



r. 

1 

-1 

1 

*2 

and “CB ■ 





• 


• 

• 


• 

X 


r 

m 


m 

. j 


_ _ 


and let N ■ 

Successively determine the rows of X as follows. First, solve the 
system Xj^B « for Xj^. Then for k - 2,3, ..., m, successively solve the 
systems 


x^B- r 


k 


k-1 

y n. . for X, . 


This concludes the algorithm. 
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The techniques utilized in the above algorithm have subsequently been 
applied by J. H. Wilkinson ("Note On The Practical Significance of The Drazln 
Inverse" Preprint) to produce a direct method for obtaining numerical solutions 
of the system of differential equations Ax(t) + Bx(t) ■ f(t) where A and B 
may both be singular matrices. 

A Perturbation Analysis Of The Limiting Probabilities For A Finite ErRodlc 
Markov Chain In Terms Of The Drazln Inverse . Given any finite homogeneous 
Markov chain with transition matrix T, the author had previously established 
the fundamental fact that the limiting probabilities were functions of the 
elements of a Drazln inverse in the following sense. If A ■ I - T, then 


11m 

n-+o 


I + T + T^ + . . . + T°~^ 
n 



where is the group Inverse of A is alternate notation for the Drazin 
inverse in the case Ind(A) 1.) The second aspect of the research conducted 
under this grant was to utilize the above relation between the group inverse 
and the limiting probabilities so as to ascertain the relative sensitivity 
the limiting probabilities might exhibit when the original chain is subject to 
a perturbation. It was discovered that a simple looking perturbation bound 
can be given in terms of A^^. The major result which was established is given 
below. 

Let T denote the transition matrix of an ergodic chain L with the limiting 
probability (row) vector w. Let L be an ergodic chain obtained from L by 
perturbing T so that T * T - E is the transition matrix of I. Let A » I - T, 

A ■ I - T, and w denote the limiting probability vector for L. The object is 


to bound the relative error 
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Although It Is not a consequence of any of the traditional theory of linear systems* 
the bound obtainable Is exactly analogous to that which appears in the analysis 
of a perturbed linear system. 

If E is "small" in the sense that || EaI^ || < 1 (for an operator norm induced 
by a vector norm (| • jj) then 



The relative error in w is bounded by the relative error in A as follows. If 

II E 11 II A^^ II < 1, then 



where G(A) ■ 1| A || [j A^'^ || . Moreover, the above Inequalities are sharp In the 
sense that there exist nontrivial examples where equality is actually attained. 

The obvious conclusion is that if G(A) is small, then the limiting probabilities 
are relatively insensitive to perturbations. If G(A) is large, the limiting 
probabilities may be quite sensitive to small perturbations. The number G(A) was 
defined by the author to be the condition of the chain . It can be demonstrated 
that G(A) also arises in dealing with the mean first passage times and other 
related aspects of a perturbed chain. 

In addition to the perturbation given above, the following perturbation 

formula was derived. If L and L are both ergodlc, then 

’ w(I + EA ) ■ w - wEA (I + EA ) , 

Independent of the size of || ||. 
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Finally I examples vere constructed to show that the bounds Involving G(A) 
give better estimates than the bounds obtainable using the t adltional theory 
based on the nonsingular "fundamental matrix" of the chain. This tends to 
reinforce the belief that the introduction of the Drazln Inverse Into the 
theory of finite Markov chains provides a significant advantage over classical 
methods. 

The results obtained under the support of this grant will be published in 
detail in four journal articles and in a portion of the author's soon to be 
released text "Generalized Inverses of Linear Transformations" (Pitman Pub. 

Co. - coauthor with S. Campbell). 



